Static and dynamic heterogeneities in a model for irreversible gelation 
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We study the structure and the dynamics in the formation of irreversible gels by means of molec- 
ular dynamics simulation of a model system where the gelation transition is due to the random 
percolation of permanent bonds between neighboring particles. We analyze the heterogeneities of 
the dynamics in terms of the fluctuations of the intermediate scattering functions: In the sol phase 
close to the percolation threshold, we find that this dynamical susceptibility increases with the time 
until it reaches a plateau. At the gelation threshold this plateau scales as a function of the wave 
vector k as k v ~ 2 , with r\ being related to the decay of the percolation pair connectedness function. 
At the lowest wave vector, approaching the gelation threshold it diverges with the same exponent 7 
as the mean cluster size. These findings suggest an alternative way of measuring critical exponents 
in a system undergoing chemical gelation. 
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In spite of the relevance of sol-gel processes in polymer 
physics, a comprehensive understanding of dynamics in 
irreversible gelation has not been achieved yet. Follow- 
ing the pioneering work by Flory and de Gennes [1], it is 
generally accepted that the divergence of the cluster size 
and the formation of a percolating network of permanent 
bonds is responsible for the sol-gel phase transition, with 
a critical increase of the viscosity coefficient and the on- 
set of an elastic response close to the gelation transition. 
The question of how the cluster size can be related to a 
thermodynamic measurable quantity like the fluctuation 
of the order parameter near a critical point, is still, to 
our knowledge, an open question. A second question of 
what are the fundamental analogies and differences be- 
tween the complex dynamics of polymer gels, colloidal 
gels, structural glasses and spin glasses has occasionally 
been discussed in the literature [2, 3], but never fully 
clarified. 

In this paper we address these two questions by means 
of theoretical arguments and molecular dynamics (MD) 
simulations. We first study by MD the gel formation in 
a model system, where neighboring particles (monomers) 
are linked by permanent bonds to form clusters of differ- 
ent sizes. By varying the volume fraction <p the system 
exhibits a percolation transition at 4> Cl in the same uni- 
versality class as random percolation. We analyze the dy- 
namics in the sol phase, by means of the self intermediate 
scattering functions and show that the percolation tran- 
sition coincides with a dynamical transition character- 
ized first by stretched exponentials and at the percolation 
threshold by a power law behavior, as found in the ex- 
periments [4] and in some recent numerical works on the 
lattice [5]. To compare chemical gelation with the slow 
dynamics observed in colloidal gels or other disordered 
systems, such as glasses and spin-glasses, we measure the 
dynamic susceptibility defined as the fluctuations of the 
self intermediate scattering function. In supercooled liq- 
uids this or similar quantities have been introduced to 



characterize the behavior of dynamical heterogeneities 
[6, 7], which typically grow with time, reach a maxi- 
mum, and then decrease at large time. This behavior, 
related to the growth of a dynamical correlation length, 
is a consequence of the transient nature of dynamic het- 
erogeneities. In this model for chemical gelation, instead, 
we find that, approaching the gelation threshold in the 
sol phase, the dynamic susceptibility increases with time, 
until it reaches a plateau in the long time limit. This 
behavior is due to the presence of static heterogeneities 
(clusters), which being persistent do not lead to the de- 
cay of the non linear susceptibility. We argue in fact 
that the non linear dynamical susceptibility, in the infi- 
nite time limit, Xas{k,(f) for k — > 0, coincides with the 
mean cluster size. The numerical data strongly support 
this result and show that, for small k, it obeys the fol- 
lowing scaling behavior Xas(k,<j>) — k v ~ 2 f(k£), where 
£ ~ (4>c — 4>)~ 1 ' is the connectedness length (the linear 
size of a critical cluster) and v the associated critical ex- 
ponent, while 2 — 77 = ~f/i/, 7 being the mean cluster size 
exponent [8] . These relations link the mean cluster size to 
the fluctuations of the intermediate scattering function, 
which can be measured from light scattering experiments 
[9]. Therefore our results can be tested experimentally 
and offer a new alternative to previous methods [10], to 
measure percolation exponents in a sol-gel transition. 

The model - We consider a 3d system of N = 1000 
particles interacting via a Lennard-Jones potential, trun- 
cated in order to have only the repulsive part: 
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where is the distance between the particles i and j. 
After a first equilibration, we introduce quenched bonds 
between particles whose relative distance is smaller than 
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Rq by adding an attractive potential: 

jjFENE = f -0.5fc -Ro H 1 - ( r ij/ R o) 2 ], Tij < R 
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representing a finitely extendable nonlinear elastic 
(FENE)[11]. The system is then further thermalized. We 
have chosen ko — 30e/<r 2 and i?o = l-5<7 as in Ref. [If] 
and performed MD simulations in a box of linear size L 
(in units of a) with periodic boundary conditions. The 
equations of motion were solved in the canonical ensem- 
ble (with a Nose-Hoover thermostat) using the velocity- 
Verlet algorithm [12] with a time step At — 0.001<5t, 
where St = cr(m/e) 1 ^ 2 , with m the mass of particle. In 
our reduced units the unit length is a, the unit energy e 
and the Boltzmann constant fee is set equal to 1. The 
temperature is fixed at T = 2 and the volume fraction 
4> = ira 3 N/6L 3 is varied from 4> = 0.02 to <j> = 0.2. 
By varying the volume fraction we find that the system 
undergoes a random percolation transition. From a stan- 
dard finite size scaling analysis [8] , we have obtained the 
percolation threshold <p c , and the critical exponents v 
(which governs the power law divergence of the connect- 
edness length £ ~ \(j> — <p c \~' J as the critical point is ap- 
proached from below) and 7 (governing the power law 
divergence of the mean cluster size x ~ \<f> ~ <^c|~ 7 )- The 
results obtained are <f) c — 0.10 ± 0.02, with critical ex- 
ponents v = 0.88 ± 0.05 and 7 — 1.8 ± 0.1 in agreement 
with random percolation. 

Dynamical properties - The dynamics at equilibrium 
is analyzed by measuring the self intermediate scattering 
function F s (k,t) defined as F s (k,t) = [(4> s (fc, t))] where 

$ a (M) = ^J2Zi eit{Pzit) ~ Pm) ^ (•■■) is the thermal 
average for a fixed bond configuration and [. . . ] is the 
average over the bond configurations. In our simulations 
the average is over 30 independent bond configurations. 
At low volume fractions F s (k,t) decays to zero follow- 
ing an exponential behavior for all the wave vectors k 
considered. Increasing the volume fraction, close to the 
percolation threshold, at low wave vectors the long time 
decay starts to follow a stretched exponential behavior 
The cluster size distribution has already 
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started to widen and therefore, over sufficiently large 
length scales (small k), the behavior of F s (k,t) is due 
to the contribution of different clusters, characterized by 
different relaxation times, whose superposition produces 
a detectable deviation from an exponential law. Close to 
4> c (see Fig. 1) the onset of a power law decay is observed 
at the lowest wave vector, k m i n = 2tt/L, indicating a crit- 
ical slowing down due to the onset of a percolating cluster 
[13]. The behavior of F s (k,t) for k = k min (plotted in 
Fig. 2 as a function of the time for different volume frac- 
tions </>) gives the relaxation dynamics over length scales 
of the order of the system size. As <f> increases towards 
4> Cl we observe a crossover from an exponential decay to 
a stretched exponential one, with decreasing as a func- 
tion of the volume fraction. At <j) c the long time decay 




FIG. 1: Fs(k,t) for decreasing values of k from left to right 
at <f> = 0.09 as a function of the time t. At k = k m i n ~ 0.35 
the data are fitted by a power law ~ t~ c with c = 0.65 ± 0.03 
(straight line). 

displays a power law behavior ~ t~ € . If (f> increases fur- 
ther, the system is out of equilibrium. For long values of 
the waiting time, we find that c decreases until the long 
time decay becomes indistinguishable from a logarithm 
behavior and eventually a two step decay appears. These 
dynamical features well reproduce the experimental ob- 
servations in different systems close to the gel transition 
[4] and agree with previous numerical results obtained on 
a lattice model [5]. 

We now analyze and discuss the behavior of 
the dynamical susceptibility X4(fc,<) associated to 
the fluctuations of $ s (fc,i), i.e. Xi{k,t) = 
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Fig. 3 for (j> < 4> c : Differently from the non monotonic 
behavior typically observed in supercooled glassy sys- 
tems, we find that X4(fc m j n ,i) increases with time until 
it reaches a plateau, whose value increases as a function 
of 0. 

Indeed we argue that in the thermodynamic limit, 
when k — > and t — > 00, %4(fc,i) tends to the mean 
cluster size. Being limt_ (00 ($ s (A; ! t)) = 0, we have in the 
long time limit 
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In Eq. (1) we have considered that, for large enough time 
t, the term e~ lfc '( r *(*) - 'j(*)) is statistically independent 
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FIG. 2: F s (k,t) for <f} = 0.02, 0.07, 0.09, 0.10, 0.11, 0.12, 0.13 
(from left to right), and k = k m i„ as a function of t. The 
full lines are, from left to right, an exponential, a stretched 
(J3 = 0.75 ±0.01 and (3 = 0.58 ±0.01) and a a power law ~ t~ c 
with c = 0.65 ± 0.03 (straight line). 

from e~* fe '( r< (°)~ r " ", so that we can factorize the ther- 
mal average. In the last term of Eq. (1), we may now 
separate the sum over connected pairs (7^ = 1, that is 
pairs belonging to the same cluster), and disconnected 
ones (jij = 0, that is pairs belonging to different clus- 
ters), so that 



lim X4(M) = ^r T 



N 



N 



7y 



-ik-(fi—rj)\ 



- T 



(1-7,,) ( 



-ik-(ri-rj) 



(2) 



For 



< 



clusters will have at most a linear size 
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nected particles will be < £. Therefore, for |fc| <C £ 1 
and jij — 1, we have (e~ lk '( Vi ~ r ^) = 1. On the other 
hand, if particles i and j are not connected, |r*i —rj\ can 
assume any value and the above argument does not ap- 
ply. However, if particles i and j are not connected, we 
can write 
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where p = N/V, /iy(f) + l = gij(r) and (1/V)gij(r) gives 
the probability density of finding the particle i in r, given 
the particle j in the origin. We have used the fact that the 
Fourier transform of 1 is zero for any wave vector different 
from zero. The correlation function hijif) decays to zero 
at a finite distance, so that (3) is of order 0(1/N). As 



FIG. 3: Main frame: X4,(kmin,t) as a function of t for 
= 0.02, 0.05, 0.06, 0.07, 0.08, 0.085, 0.09, 0.095, 0.10 (from 
bottom to top). Inset: Asymptotic values of the susceptibil- 
ity (full triangles), Xas(k m in, 4>) an d mean cluster size (open 
squares) as a function of {4> c — (f>). The data are fitted by the 
power law (<f) c — 0)~ 7 with 7 = 1.8 ± 0.1. 



there are at most N 2 disconnected pairs, the second term 
in Eq. (2) is of order 0(l/N), so that it can be neglected 
in the thermodynamic limit. Finally, we have 
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In the inset of Fig. 3 the asymptotic value Xas{k m in, 4>) = 
linit^oo Xii^min, t) is plotted as a function of {4> c — 4>) 
together with the mean cluster size. We find that, 
as the percolation threshold is approached from below, 
Xas(k m in,4>) shows a power law behavior with an expo- 
nent which, within the numerical accuracy, is in agree- 
ment with the value of the exponent 7 of the mean cluster 
size. 

Using a standard scaling argument [8], we can easily 
predict the behavior of Xas(k,(f)) for small fc close to (j) c . 
Since k scales as the inverse of a length, Xas(k,cf)) = 
k r '- 2 f{k^), where £ ~ (4> c - 4>)- v , 2 - 7/ = 7/1/ and 
f(z) ~ z 1 ^ for small values z and f(z) ~ const for 
large values of z. Fig. 4 shows that the numerical data 
strongly support this scaling behavior. 

In conclusion according to these results the percolation 
critical exponents in the sol-gel transition could be mea- 
sured via the fluctuations of the self intermediate scat- 
tering functions by means of light scattering. Moreover 
our results confirm that one key difference between irre- 
versible gelation due to chemical bonds and supercooled 
liquids close to the glass transition is that in irreversible 
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FIG. 4: Main frame: Xas(k,4>c) as a function of the wave 
vector k and, from bottom to top, <j> = 0.02, 0.05, 0.06, 0.07, 
0.08, 0.085, 0.09, 0.095. At cj> c the data are fitted with ~ 
(fc)"" 2 (full line), with 2 - r\ = 2.03 ± 0.02. Inset: Scaling 
plot of Xas(k, 4>)k 2 ~ v as a function of fc _1 ^(<^ c — <f>). 

gelation the heterogeneities have a static nature (clus- 
ters). These clusters, on the other hand, affect the dy- 
namics and as a consequence the dynamic transition co- 
incides with the static transition, characterized by the di- 
vergence of a static correlation length (linear size of the 
clusters). In this respect the dynamical slowing down 
here is similar to the critical slowing down found close 
to a second order critical point, with the relaxation time 
diverging with a power law in the control parameters. In- 
terestingly the behavior of the dynamical susceptibility 
is very similar to the one observed in spin glass models 
[14], where for long times it tends to the non linear static 
susceptibility, which diverges at the spin glass critical 
point. In fact, also in this case the dynamic transition is 
connected to the divergence of a static length. 

What can we expect in colloidal gelation, where the 
bonds are not quenched and have a finite lifetime? In 
principle, due to finite bond lifetime, the clusters are not 
permanent anymore, consequently we expect a structural 
arrest with a dynamical susceptibility of the type found 
in glass forming liquid. This is in fact found in experi- 
mental investigations of colloidal suspension [15] and in 
some molecular dynamic simulations [16]. Interestingly 
enough, also in some spin glass models, by introducing 
interactions with finite lifetime, there is no divergence of 
the static susceptibility and one recovers the behavior of 
the dynamic susceptibility which is typical of supercooled 
liquids [14]. 
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